rm(list=ls())

library(tidyverse)
library(BBmisc)
library(data.table)
library(rmeta)
source('Habit_source.R', echo=TRUE)


l = rbindlist(lapply(seq(14, 352,  by = 14), FUN = function(kt){
output = list.files('results')[grepl(paste0("bw",kt,".Rdata"), list.files('results'))  ]

main_results = as.data.frame(matrix(nrow = 0, ncol =6))
colnames(main_results) = c("coef", "se", "race", "dv",  "endog", "state")



for(check in output){
  load(paste0("results/", check))
  results = as.data.frame(matrix(nrow = 21, ncol = 6))
  colnames(results) = c("coef", "se", "race", "dv", "endog", "state")
  
  ## recode if error
  
  # white
  if(is.error(ivmod1_white)){
    ivmod1_white = list()
    ivmod1_white$cace=NA
    ivmod1_white$se=NA
  }
  
  if(is.error(ivmod2_white)){
    ivmod2_white = list()
    ivmod2_white$cace=NA
    ivmod2_white$se=NA
  }
  
  if(is.error(ivmod5_white)){
    ivmod5_white = list()
    ivmod5_white$cace=NA
    ivmod5_white$se=NA
  }
  
  if(is.error(ivmod6_white)){
    ivmod6_white = list()
    ivmod6_white$cace=NA
    ivmod6_white$se=NA
  }
  
  if(is.error(ivmod_rep_white)){
    ivmod_rep_white = list()
    ivmod_rep_white$cace=NA
    ivmod_rep_white$se=NA
  }
  

  #black
  if(is.error(ivmod1_black)){
    ivmod1_black = list()
    ivmod1_black$cace=NA
    ivmod1_black$se=NA
  }
  
  if(is.error(ivmod2_black)){
    ivmod2_black = list()
    ivmod2_black$cace=NA
    ivmod2_black$se=NA
  }
  
  if(is.error(ivmod5_black)){
    ivmod5_black = list()
    ivmod5_black$cace=NA
    ivmod5_black$se=NA
  }
  
  if(is.error(ivmod6_black)){
    ivmod6_black = list()
    ivmod6_black$cace=NA
    ivmod6_black$se=NA
  }
  
  if(is.error(ivmod_rep_black)){
    ivmod_rep_black = list()
    ivmod_rep_black$cace=NA
    ivmod_rep_black$se=NA
  }
  

  
  #hispanic
  if(is.error(ivmod1_hispanic)){
    ivmod1_hispanic = list()
    ivmod1_hispanic$cace=NA
    ivmod1_hispanic$se=NA
  }
  
  if(is.error(ivmod2_hispanic)){
    ivmod2_hispanic = list()
    ivmod2_hispanic$cace=NA
    ivmod2_hispanic$se=NA
  }
  
  if(is.error(ivmod5_hispanic)){
    ivmod5_hispanic = list()
    ivmod5_hispanic$cace=NA
    ivmod5_hispanic$se=NA
  }
  
  if(is.error(ivmod6_hispanic)){
    ivmod6_hispanic = list()
    ivmod6_hispanic$cace=NA
    ivmod6_hispanic$se=NA
  }
  
  if(is.error(ivmod_rep_hispanic)){
    ivmod_rep_hispanic = list()
    ivmod_rep_hispanic$cace=NA
    ivmod_rep_hispanic$se=NA
  }
  

  ##
  results$coef = c(ivmod1_white$cace, ivmod1_black$cace, ivmod1_hispanic$cace,
                   ivmod2_white$cace, ivmod2_black$cace, ivmod2_hispanic$cace,
                   ivmod3_white$cace, ivmod3_black$cace, ivmod3_hispanic$cace,
                   ivmod4_white$cace, ivmod4_black$cace, ivmod4_hispanic$cace,
                   
                   ivmod5_white$cace, ivmod5_black$cace, ivmod5_hispanic$cace,
                   ivmod6_white$cace, ivmod6_black$cace, ivmod6_hispanic$cace,
                   ivmod_rep_white$cace, ivmod_rep_black$cace, ivmod_rep_hispanic$cace)
  
  results$se = c(  ivmod1_white$se, ivmod1_black$se, ivmod1_hispanic$se,
                   ivmod2_white$se, ivmod2_black$se, ivmod2_hispanic$se,
                   ivmod3_white$se, ivmod3_black$se, ivmod3_hispanic$se,
                   ivmod4_white$se, ivmod4_black$se, ivmod4_hispanic$se,
                   
                    ivmod5_white$se, ivmod5_black$se, ivmod5_hispanic$se,
                   ivmod6_white$se, ivmod6_black$se, ivmod6_hispanic$se,
                   ivmod_rep_white$se, ivmod_rep_black$se, ivmod_rep_hispanic$se)
  

   ##############
  results$race = c(rep(c("white", "black", "hispanic"), 7))
  
  results$dv = c(c(rep("Vote 2016", 3),
                 rep("Vote 2014", 3),
                 rep("Vote 2012", 3),
                 rep("Vote 2010", 3),
                 rep("Democrat", 3),
                 rep("Registered", 3),
                 rep("Republican", 3)))
  
  
  results$endog = c(rep("Vote 2008", 3*7))
  results$state = sub(as.character(kt),"",sub("_bw","",sub("results_","", sub(".Rdata", "", check))))
  results$bw = kt
  
  main_results = rbind(main_results, results)
  rm(list=ls(pattern = "mod"))
  
}





print(kt)
return(main_results)
}))


m = rbindlist(lapply(seq(14, 352,  by = 14),FUN=function(kt){
  data = l %>%
    filter(bw == kt & !is.na(coef) & !is.na(se))
predictions_white = as.data.frame(matrix(NA, ncol = 5))
predictions_black = as.data.frame(matrix(NA, ncol = 5))
predictions_hispanic = as.data.frame(matrix(NA, ncol = 5))
colnames(predictions_white) = c("ev", "dv", "coef", "se", "race")
colnames(predictions_black) = c("ev", "dv", "coef", "se", "race")
colnames(predictions_hispanic) = c("ev", "dv", "coef", "se", "race")
predictions_white$race = "white"
predictions_black$race = "black"
predictions_hispanic$race = "hispanic"

for(i in 1:length(unique(data$endog))){
  for(j in 1:length(unique(data$dv))){
    if(length(data$dv[data$dv == unique(data$dv)[j] & data$endog==unique(data$endog)[i]]) > 0){
  
      preds_white = meta.maker(data$coef[data$endog == unique(data$endog)[i] & data$race == "white" & data$dv == unique(data$dv)[j] ], data$se[data$endog == unique(data$endog)[i] & data$race == "white"& data$dv == unique(data$dv)[j]])
      predictions_white = rbind(predictions_white, c(unique(data$endog)[i], unique(data$dv)[j], preds_white$est, preds_white$se, "white"))
      
      preds_black = meta.maker(data$coef[data$endog == unique(data$endog)[i] & data$race == "black" & data$dv == unique(data$dv)[j] ], data$se[data$endog == unique(data$endog)[i] & data$race == "black"& data$dv == unique(data$dv)[j]])
      predictions_black = rbind(predictions_black, c(unique(data$endog)[i], unique(data$dv)[j], preds_black$est, preds_black$se, "black"))
      
      preds_hispanic = meta.maker(data$coef[data$endog == unique(data$endog)[i] & data$race == "hispanic" & data$dv == unique(data$dv)[j]], data$se[data$endog == unique(data$endog)[i] & data$race == "hispanic"& data$dv == unique(data$dv)[j]])
      predictions_hispanic = rbind(predictions_hispanic, c(unique(data$endog)[i], unique(data$dv)[j], preds_hispanic$est, preds_hispanic$se, "hispanic"))
    }
    
  }
  
}
meta_plot=rbind(predictions_white, predictions_black, predictions_hispanic)
for(i in 3:(ncol(meta_plot)-1)){
  meta_plot[,i] = as.numeric(meta_plot[,i])
}
meta_plot$bw=kt
return(meta_plot)
}))



ggsave(plot=ggplot(data = m %>%
                     filter(grepl('Vote',dv)), aes(x=bw, y=coef, group = race), color = "black") +
         geom_line(aes(linetype=race))+
         geom_hline(yintercept = 0, linetype = 'dotted', color = 'grey60') +
         geom_ribbon(aes(ymin = coef-1.96*se, ymax=coef+1.96*se), alpha = .25)+
         xlab("Bandwidth")+
         ylab("CACE")+
         theme_bw()+
         theme(text = element_text(size=15), panel.grid.major = element_blank(), panel.grid.minor = element_blank(),legend.position="bottom", 
               legend.title = element_blank(),
               legend.text=element_text(size=15), axis.text = element_text(size = 14),
               legend.key.width=unit(3,"line")) +
         facet_wrap(~dv), 
       filename = "figures/FigS11.pdf", dpi = 300, width = 6, height = 6, units = 'in')


ggsave(plot=ggplot(data = m %>%
                     filter(!grepl('Vote',dv) & !is.na(dv)), aes(x=bw, y=coef, group = race), color = "black") +
         geom_line(aes(linetype=race))+
         geom_hline(yintercept = 0, linetype = 'dotted', color = 'grey60') +
         geom_ribbon(aes(ymin = coef-1.96*se, ymax=coef+1.96*se), alpha = .25)+
         xlab("Bandwidth")+
         ylab("CACE")+
         theme_bw()+
         theme(text = element_text(size=15), panel.grid.major = element_blank(), panel.grid.minor = element_blank(),legend.position="bottom", 
               legend.title = element_blank(),
               legend.text=element_text(size=15), axis.text = element_text(size = 14),
               legend.key.width=unit(3,"line")) +
         facet_wrap(~dv), 
       filename = "figures/FigS12.pdf", dpi = 300, width = 6, height = 3, units = 'in')





